---
title: "APSR Replication"
author: "Velez, Liu"
output:
  html_document: default
  pdf_document: default
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(warning = FALSE, echo = FALSE, message = FALSE)
```

# Main Figures

## Figure 1
```{r}
source("scripts/clean_study1.R")

study1 %>%
 group_by(motivation) %>%
 do(mod = broom::tidy(lm_robust(post_attitudes ~ pre_attitudes + 
                                 information, .,
                                clusters = id))) %>%
 unnest(c(mod)) %>%
 filter(term != "pre_attitudes" & term != "(Intercept)") %>%
 mutate(term = plyr::mapvalues(term, c("informationcon", "informationmixed"),
                               c("Con", "Mixed")),
        Motivation = plyr::mapvalues(motivation,
                                     c("accuracy", "directional"),
                                     c("Accuracy Prime",
                                       "Directional Prime"))) %>%
 mutate(term = forcats::fct_rev(term)) %>%
 ggplot(aes(x = term, y = estimate, shape = Motivation)) +
 geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                 position = position_dodge(.25),
                 size = 1.25) +
 geom_linerange(aes(ymin = estimate - 1.405*std.error,
                    ymax = estimate + 1.405*std.error),
                position = position_dodge(.25),
                size = 2.5) +
 annotate("text", x = 0.5, y = .125, label = "Polarization",
          angle = 90, size = 5, family = "Roboto Condensed") +
 annotate("text", x = 0.5, y = -.125, label = "Moderation", 
          angle = 90, size = 5, family = "Roboto Condensed") +
 annotate("text", x = 0.5, y = .23, label = "\u2192",
          angle = 90, size = 5, family = "Arial") +
 annotate("text", x = 0.5, y = -.23, label = "\u2190", 
          angle = 90, size = 5, family = "Arial") +
 geom_hline(yintercept = 0, linetype = "dashed") +
 theme_ipsum_rc() +
 labs(y = "Attitude Strength (1-7)", x = "") +
 theme_params +
 theme(legend.position = "none") +
 ylim(-.25, .25) 

study1 %>%
 group_by(motivation) %>%
 do(mod = broom::tidy(lm_robust(post_certainty ~ pre_certainty + information, .,
                                clusters = id))) %>%
 unnest(c(mod)) %>%
 filter(term != "pre_certainty" & term != "(Intercept)") %>%
 mutate(term = plyr::mapvalues(term, c("informationcon", "informationmixed"),
                               c("Con", "Mixed")),
        Motivation = plyr::mapvalues(motivation,
                                     c("accuracy", "directional"),
                                     c("Accuracy Prime",
                                       "Directional Prime"))) %>%
 mutate(term = forcats::fct_rev(term)) %>%
 ggplot(aes(x = term, y = estimate, shape = Motivation)) +
 geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                 position = position_dodge(.25),
                 size = 1.25) +
 geom_linerange(aes(ymin = estimate - 1.405*std.error,
                     ymax = estimate + 1.405*std.error),
                 position = position_dodge(.25),
                 size = 2.5) +
 geom_hline(yintercept = 0, linetype = "dashed") +
 theme_ipsum_rc() +
 labs(y = "Attitude Certainty (0-100)", x = "") +
 theme_params +
 ylim(-2, 2) 


```

## Figure 2
```{r}
source("scripts/clean_study2.R")

combined_df %>%
  filter(`Issue Strength` == "Strong") %>%
  group_by(motivation) %>%
  do(mod = broom::tidy(emmeans(lm_robust(post_attitudes ~ information +
                                           pre_attitudes, .,
                                         clusters = participantId),
                               specs = pairwise ~ information)$emmeans)) %>%
  unnest(c(mod)) %>%
  ggplot(aes(x = information, y = estimate,
             group = motivation)) +
  geom_point(position = position_dodge(.5)) +
  geom_point(position = position_dodge(.5), size = 4) +
  geom_linerange(aes(ymin = estimate - 1.405*std.error,
                     ymax = estimate + 1.405*std.error),
                 position = position_dodge(.5), size = 2) +
  geom_linerange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                 position = position_dodge(.5)) +
  theme_ipsum_rc() +
  facet_wrap(. ~ motivation, scales = "free") +
  labs(y = "Estimate") +
  theme(panel.background = element_rect(fill = 'gray95'),
        panel.grid.minor = element_blank(),
        panel.grid.major = element_blank(),
        axis.text.x = element_text(size = 15),
        axis.text.y = element_text(size = 15),
        axis.title.y = element_text(size = 15),
        axis.title.x = element_text(size = 15),
        legend.position = "bottom") +
  labs(x = "", y = "Attitude Strength (1-7)") +
  scale_color_grey(start = 0, end = .6)

combined_df %>%
  filter(`Issue Strength` == "Strong") %>%
  group_by(motivation) %>%
  do(mod = broom::tidy(emmeans(lm_robust(post_certain_scores ~ information +
                                           pre_certain_scores, .,
                                         clusters = participantId),
                               specs = pairwise ~ information)$emmeans)) %>%
  unnest(c(mod)) %>%
  mutate(information = plyr::mapvalues(information,
                                       c("con", "mixed", "pro"),
                                       c("Cons",
                                         "Mixed",
                                         "Pros")),
         motivation = plyr::mapvalues(motivation,
                                      c("accuracy", "directional"),
                                      c("Accuracy Prime",
                                        "Directional Prime"))) %>%
  ggplot(aes(x = information, y = estimate)) +
  geom_point(position = position_dodge(.5)) +
  geom_point(position = position_dodge(.5), size = 4) +
  geom_linerange(aes(ymin = estimate - 1.405*std.error,
                     ymax = estimate + 1.405*std.error),
                 position = position_dodge(.5), size = 2) +
  geom_linerange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                 position = position_dodge(.5)) +
  theme_ipsum_rc() +
  facet_wrap(. ~ motivation,
             scales = "free", ncol = 4) +
  labs(y = "Estimate") +
  theme(panel.background = element_rect(fill = 'gray95'),
        panel.grid.minor = element_blank(),
        panel.grid.major = element_blank(),
        axis.text.x = element_text(size = 15),
        axis.text.y = element_text(size = 15),
        axis.title.y = element_text(size = 15),
        axis.title.x = element_text(size = 15),
        legend.position = "bottom") +
  labs(x = "", y = "Attitude Certainty (1-9)") +
  scale_color_grey(start = 0, end = .6)
```

## Figure 3
```{r}
source("scripts/clean_study3.R")

study3 %>%
  lm_robust(post_outcome_1 ~ treatment + age + education + hhi
            + ideology + white + political_party,
            .) %>%
  broom::tidy() %>%
  filter(term == "treatment") %>%
  ggplot(aes(x = term, y = estimate)) +
  geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                      ymax = estimate + 1.96*std.error),
                  size = 1.25) +
  ylim(-.70, .70) +
  theme_ipsum_rc() +
  theme_params +
  theme(legend.position = "none",
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        axis.text.x = element_blank()) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  labs(x = "Attitude Strength (1-7)", y = "Average Treatment Effect") +
  scale_x_discrete(position = "top")

study3 %>%
  lm_robust(post_certain_scores ~ treatment + age + education + hhi
            + ideology + white + political_party,
            .) %>%
  broom::tidy() %>%
  filter(term == "treatment") %>%
  ggplot(aes(x = term, y = estimate)) +
  geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                      ymax = estimate + 1.96*std.error),
                  size = 1.25) +
  ylim(-.7, .7) +
  theme_ipsum_rc() +
  theme_params +
  theme(legend.position = "none",
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        axis.text.x = element_blank())  +
  labs(x = "Attitude Certainty (1-9)", y = "") +
  geom_hline(yintercept = 0, linetype = "dashed") +
  scale_x_discrete(position = "top")

```

## Figure 4
```{r}
source("scripts/clean_study4.R")

study4_coded %>%
 # First model
 do(mod = broom::tidy(lm_robust(pw2p_scale ~ treatment + likert + certainty + age + educ + income + ideology +
                              participation_scale + pk_scale + pid7,
                             data = ., subset = is.na(generic_flag)))) %>%
 bind_rows(study4_coded %>%
            do(mod = broom::tidy(lm_robust(pw2p_defend ~ treatment + likert + certainty + age + educ + income + ideology +
                                         participation_scale + pk_scale + pid7,
                                        data = ., subset = is.na(generic_flag))))) %>%
 mutate(outcome_label = c("Attitude Extremity",
                          "Attitude Defense")) %>%
 # Filter terms
 unnest(c(mod)) %>%
 filter(term == "treatment") %>%
 # Plotting
 ggplot(aes(x = term,
            y = estimate)) +
 geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                  position = position_dodge(.5), size = 1.2) +
 geom_hline(aes(yintercept = 0), linetype = "dashed") +
 theme_ipsum_rc() +
 labs(y = "Average Treatment Effect", x = "") +
 theme_params +
 facet_grid(. ~ outcome_label) +
 theme(axis.text.x = element_blank()) +
 ylim(c(-.35, .35)) +
 annotate("text", x = 0.5, y = .125, label = "Polarization",
          angle = 90, size = 5, family = "Roboto Condensed") +
 annotate("text", x = 0.5, y = -.125, label = "Moderation",
          angle = 90, size = 5, family = "Roboto Condensed") +
 annotate("text", x = 0.5, y = .23, label = "\u2192",
          angle = 90, size = 5, family = "Arial") +
 annotate("text", x = 0.5, y = -.23, label = "\u2190",
          angle = 90, size = 5, family = "Arial")
```

## Figure 5
```{r}
source("scripts/clean_study5.R")

study5_coded %>%
  # First model
  do(mod = broom::tidy(lm_robust(pw2p_scale ~ as.factor(treatment),
                                 + likert + certainty + age + educ + income + ideology +
                                   participation_scale + pk_scale + pid7,
                                 data = ., subset = is.na(generic_flag)))) %>%
  # Second model
  bind_rows(study5_coded %>%
              do(mod = broom::tidy(lm_robust(pw2p_defend ~ as.factor(treatment),
                                             + likert + certainty + age + educ + income + ideology +
                                               participation_scale + pk_scale + pid7,
                                             data = ., subset = is.na(generic_flag))))) %>%
  # Adding outcome labels
  mutate(outcome_label = c("Attitude Extremity Scale",
                           "Attitude Defense Scale")) %>%
  # Filter terms
  unnest(c(mod)) %>%
  filter(grepl("treatment", term)) %>%
  mutate(term = plyr::mapvalues(term,
                                c("as.factor(treatment)1", "as.factor(treatment)2"),
                                c("Low Valence\nCounterargument",
                                  "High Valence\nCounterargument"))) %>%
  # Plotting
  ggplot(aes(x = forcats::fct_reorder(term,estimate),
             y = estimate,
             shape = term)) +
  geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                      ymax = estimate + 1.96*std.error),
                  position = position_dodge(.5), size = 1) +
  geom_linerange(aes(ymin = estimate - 1.405*std.error,
                     ymax = estimate + 1.405*std.error),
                 position = position_dodge(.5), size = 3) +
  geom_hline(aes(yintercept = 0), linetype = "dashed") +
  theme_ipsum_rc() +
  facet_grid(. ~ outcome_label, scales = "free") +
  labs(y = "Average Treatment Effect", x = "") +
  theme(panel.background = element_rect(fill = 'gray95'),
        panel.grid.minor = element_blank(),
        panel.grid.major = element_blank(),
        legend.text = element_text(size = 15),
        axis.text.x = element_text(size = 15),
        legend.title=element_blank(),
        axis.text.y = element_text(size = 15),
        axis.title.y = element_text(size = 15),
        axis.title.x = element_text(size = 15),
        legend.position = "none",
        plot.margin = unit(c(.25,.25,.25,.25), "cm"),
        strip.text.x = element_text(size = 15))
```

# Appendix

## Appendix A.1 (Validation Tests)

```{r}
source("scripts/clean_study1.R")

library(modelsummary)

batch <- rio::import("data/mturk_codes.csv")

study1 %>%
 left_join(batch %>% 
             transmute(id = Input.id,
                       denigrating_comments = as.numeric(plyr::mapvalues(Answer.intent.label,
                                                                         c("All responses agree with the prompt",
                                                                           "1 response rejects/disagrees with the prompt",
                                                                           "2 responses reject/disagree with the prompt",
                                                                           "3 responses reject/disagree with the prompt",
                                                                           "4 responses reject/disagree with the prompt"),
                                                                         c(0,1,2,3,4)))) %>% 
             group_by(id) %>%
             summarize(denigrating_comments = mean(denigrating_comments,na.rm=T))) %>%
 left_join(mturk %>% transmute(id, n), by = "id") %>%
 mutate(share_den = denigrating_comments/n) -> test_r

modelsummary(list(lm_robust(share_den ~ information, test_r),
     lm_robust(share_den ~ information * motivation, test_r),
     lm_robust(log(timer) ~ information, test_r),
     lm_robust(log(timer) ~ information * motivation, test_r)))

```

## Appendix A.2 (Differences between Strong and Weak Attitudes)

```{r}
source("scripts/clean_study2.R")

combined_df %>%
 group_by(motivation, `Issue Strength`) %>%
 do(mod = broom::tidy(lm_robust(post_attitudes ~ pre_attitudes + information, .,
                                clusters = participantId))) %>%
 unnest(c(mod)) %>%
 filter(term != "pre_attitudes" & term != "(Intercept)") %>%
 mutate(term = plyr::mapvalues(term, c("informationCon", "informationMixed"),
                               c("Con", "Mixed")),
        Motivation = plyr::mapvalues(motivation,
                                     c("accuracy", "directional"),
                                     c("Accuracy Prime",
                                       "Directional Prime"))) %>%
 mutate(term = forcats::fct_rev(term)) %>%
 ggplot(aes(x = term, y = estimate, shape = Motivation)) +
 geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                 size = 1.25, position = position_dodge(.5)) +
 geom_linerange(aes(ymin = estimate - 1.405*std.error,
                    ymax = estimate + 1.405*std.error),
                position = position_dodge(.5),
                size = 2.5) +
 facet_wrap(. ~ `Issue Strength`) +
 geom_hline(yintercept = 0, linetype = "dashed") +
 theme_ipsum_rc() +
 labs(y = "Attitude Strength (1-7)", x = "") +
 theme_params +
 ylim(-.8, .8)

combined_df %>%
 group_by(motivation, `Issue Strength`) %>%
 do(mod = broom::tidy(lm_robust(post_certain_scores ~ 
                                 pre_certain_scores + information, .,
                                clusters = participantId))) %>%
 unnest(c(mod)) %>%
 filter(term != "pre_certain_scores" & term != "(Intercept)") %>%
 mutate(term = plyr::mapvalues(term, c("informationCon", "informationMixed"),
                               c("Con", "Mixed")),
        Motivation = plyr::mapvalues(motivation,
                                     c("accuracy", "directional"),
                                     c("Accuracy Prime",
                                       "Directional Prime"))) %>%
 mutate(term = forcats::fct_rev(term)) %>%
 ggplot(aes(x = term, y = estimate, shape = Motivation)) +
 geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                 size = 1.25, position = position_dodge(.5)) +
 geom_linerange(aes(ymin = estimate - 1.405*std.error,
                    ymax = estimate + 1.405*std.error),
                position = position_dodge(.5),
                size = 2.5) +
 facet_wrap(. ~ `Issue Strength`, scales = "free") +
 geom_hline(yintercept = 0, linetype = "dashed") +
 theme_ipsum_rc() +
 theme(panel.grid.major = element_blank()) +
 labs(y = "Attitude Certainty (1-9)", x = "")  +
 theme_params +
 ylim(-1, 1) 

```

## Appendix A.3 (Prior Attitude Effect)

```{r}
combined_df %>%
 group_by(condition) %>%
 do(mod = broom::tidy(emmeans(lm_robust(post_con_strength ~ motivation + 
                                         pre_attitudes, .,
                                        clusters = participantId),
                              specs = pairwise ~ motivation)$emmeans)) %>%
 unnest(c(mod)) %>%
 bind_cols(outcome = "Strength Rating",
           type = "Con Arguments") %>%
 bind_rows(combined_df %>%
            group_by(condition) %>%
            do(mod = broom::tidy(emmeans(lm_robust(post_pro_strength ~ motivation + 
                                                    pre_attitudes, .,
                                                   clusters = participantId),
                                         specs = pairwise ~ motivation)$emmeans)) %>%
            unnest(c(mod)) %>%
            bind_cols(outcome = "Strength Rating",
                      type = "Pro Arguments")) %>%
 bind_rows(combined_df %>%
            group_by(condition) %>%
            do(mod = broom::tidy(emmeans(lm_robust(post_con_acc ~ motivation + 
                                                    pre_attitudes, .,
                                                   clusters = participantId),
                                         specs = pairwise ~ motivation)$emmeans)) %>%
            unnest(c(mod)) %>%
            bind_cols(outcome = "Accuracy Rating",
                      type = "Con Arguments")) %>%
 bind_rows(combined_df %>%
            group_by(condition) %>%
            do(mod = broom::tidy(emmeans(lm_robust(post_pro_acc ~ motivation + 
                                                    pre_attitudes, .,
                                                   clusters = participantId),
                                         specs = pairwise ~ motivation)$emmeans)) %>%
            unnest(c(mod)) %>%
            bind_cols(outcome = "Accuracy Rating",
                      type = "Pro Arguments")) %>%
 mutate(`Issue Strength` = plyr::mapvalues(condition,
                                           c("strong", "weak"),
                                           c("Strong",
                                             "Weak")),
        motivation = plyr::mapvalues(motivation,
                                     c("accuracy", "directional"),
                                     c("Accuracy Prime",
                                       "Directional Prime"))) %>%
 ggplot(aes(x = `Issue Strength`, y = estimate, 
            color = motivation,
            shape = type)) +
 geom_point(position = position_dodge(.5), size = 6) +
 geom_linerange(aes(ymin = estimate - 1.405*std.error,
                    ymax = estimate + 1.405*std.error),
                position = position_dodge(.5), size = 1.25) +
 geom_linerange(aes(ymin = estimate - 1.96*std.error,
                    ymax = estimate + 1.96*std.error),
                position = position_dodge(.5), size = .25) +
 theme_ipsum_rc() +
 facet_wrap(. ~ outcome, scales = "free", nrow = 2) +
 coord_flip()+
 labs(y = "Estimate") +
 theme_params +
 labs(x = "", y = "Ratings") +
 scale_color_grey(start = 0, end = .6)
```

## Appendix A.4 (Disconfirmation Bias)

```{r}
source("scripts/clean_study2.R")

combined_df %>%
 group_by(`Issue Strength`) %>%
 do(mod = broom::tidy(lm_robust(share_den ~ information, .))) %>%
 unnest(c(mod)) %>%
 mutate(term = plyr::mapvalues(term, c("informationCon", "informationMixed"),
                               c("Con", "Mixed"))) %>%
 filter(!grepl("Intercept|pre", term)) %>%
 ggplot(aes(x = forcats::fct_rev(term), y = estimate,
            color = `Issue Strength`)) +
 geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                 position = position_dodge(.5),
                 size = 1.2) +
 geom_linerange(aes(ymin = estimate - 1.405*std.error,
                     ymax = estimate + 1.405*std.error),
                 position = position_dodge(.5),
                 size = 2) +
 theme_ipsum_rc() +
 theme_params +
 scale_color_grey(start = 0, end = .6) +
 geom_hline(yintercept = 0, linetype = "dashed") +
 labs(x = "", y = "Effect on Denigrating Share") + 
 theme(legend.position = "none") 

combined_df %>%
 group_by(`Issue Strength`) %>%
 do(mod = broom::tidy(lm_robust(log(timer) ~ information, .))) %>%
 unnest(c(mod)) %>%
 mutate(term = plyr::mapvalues(term, c("informationCon", "informationMixed"),
                               c("Con", "Mixed"))) %>%
 filter(!grepl("Intercept|pre", term)) %>%
 ggplot(aes(x = forcats::fct_rev(term), y = estimate,
            color = `Issue Strength`)) +
 geom_pointrange(aes(ymin = estimate - 1.96*std.error,
                     ymax = estimate + 1.96*std.error),
                 position = position_dodge(.5),
                 size = 1.2) +
 geom_linerange(aes(ymin = estimate - 1.405*std.error,
                    ymax = estimate + 1.405*std.error),
                position = position_dodge(.5),
                size = 2) +
 theme_ipsum_rc() + 
 theme_params +
 geom_hline(yintercept = 0, linetype = "dashed") +
 labs(x = "", y = "Effect on Time (Logged)") +
 scale_color_grey(start = 0, end = .6) +
 theme(legend.position = "bottom") 
```

## Appendix B.1 (Distribution of Issue Subtopics)

```{r}
source("scripts/prep_appendix_B1.R")

stats_subtopic %>% 
 ggplot(aes(x = prop, y = reorder(subtopic, prop))) + 
 geom_point(aes(shape = study, color = study), size = 2.5) +
 facet_wrap(~ chart, nrow = 1) +
 hrbrthemes::theme_ipsum_rc() +
 theme(panel.background = element_rect(fill = 'gray95'),
       panel.grid.minor = element_blank(),
       panel.grid.major = element_blank(),
       axis.text.x = element_text(size = 10),
       axis.text.y = element_text(size = 8.7),
       axis.title.y = element_text(size = 11),
       axis.title.x = element_text(size = 11),
       legend.position = "bottom") +
 labs(x = "Proportion of Responses", y = "Participant-Generated Subtopics") +
 scale_color_manual(name="", 
                    labels = c("Experiment 1", "Experiment 2", "Experiment 3", "Total"),
                    values = c("royalblue2", "olivedrab3", "tomato", "gray20")) +
 scale_shape_manual(name = "",
                    labels = c("Experiment 1", "Experiment 2", "Experiment 3", "Total"),
                    values = c(21, 24, 22, 18))
```

## Appendix B.2 (Are GPT-3 Arguments Persuasive?)

```{r}
source("scripts/prep_appendix_B2.R")
library(estimatr)

m1 <- lm_robust(rating ~ 1, combined_pairs, clusters = participant)
m2 <- lm_robust(rating ~ valence, combined_pairs, clusters = participant)

list(m1, m2) %>%
  modelsummary()

combined_pairs %>%
  group_by(subtopic) %>%
  do(mod = broom::tidy(lm_robust(rating ~ 1, .,
                                 clusters = participant))) %>%
  unnest(c(mod)) %>%
  ggplot(aes(y = forcats::fct_reorder(subtopic, estimate),
             x = estimate)) +
  geom_pointrange(aes(xmin = estimate - 1.96*std.error,
                      xmax = estimate + 1.96*std.error),
                  position = position_dodge(.5),
                  size = 1.10) +
  theme_ipsum_rc() +
  theme(panel.grid.major = element_blank()) +
  geom_vline(xintercept = .5, linetype = "dashed") +
  labs(x = "Estimate", y = "")

combined_pairs %>%
  group_by(subtopic) %>%
  do(mod = broom::tidy(lm_robust(rating ~ 1, .,
                                 clusters = participant))) %>%
  unnest(c(mod)) %>% 
  select(subtopic, estimate, std.error, p.value, conf.low, conf.high) %>% 
  mutate(estimate = round(estimate, 3),
         std.error = round(std.error, 3),
         p.value = round(p.value, 4),
         conf.low = round(conf.low, 3),
         conf.high = round(conf.high, 3)) %>% 
  xtable::xtable()
```

## Appendix B.3 (Visualizing Ceiling Effects)

```{r}
source("scripts/prep_appendix_B3.R")

# Figure B8
pre_control %>% 
  filter(study!="Study 3" | (study=="Study 3" & treatment==0)) %>% 
  filter(!is.na(strength)) %>% 
  group_by(study, n, strength) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 3", 959, n), # adjust n in S3 to control group
         prop = freq/n) %>%
  ggplot(aes(x = factor(strength), y = prop)) + 
  geom_bar(stat = "identity", width = 0.7, alpha = 0.7) + 
  facet_grid(study ~ .) +
  theme_params + 
  labs(x = "Attitude Strength (1-7)", y = "Proportion") +
  theme(strip.text = element_text(size = 7),
        axis.text.y = element_text(size = 7),
        axis.title.y = element_text(size = 7),
        axis.text.x = element_text(size = 7),
        axis.title.x = element_text(size = 7))

pre_control %>% 
  filter(study!="Study 3" | (study=="Study 3" & treatment==0)) %>% 
  filter(!is.na(certainty)) %>% 
  group_by(study, n, certainty) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 3", 916, n), # adjust n in S3 to control group
         prop = freq/n) %>%
  ggplot(aes(x = certainty, y = prop)) + 
  geom_bar(stat = "identity", width=0.85, alpha = 0.7) + 
  facet_grid(study ~ .) +
  theme_params + 
  labs(x = "Attitude Certainty (0-100)", y = "Proportion") +
  theme(strip.text = element_text(size = 7),
        axis.text.y = element_text(size = 7),
        axis.title.y = element_text(size = 7),
        axis.text.x = element_text(size = 7),
        axis.title.x = element_text(size = 7))

# Figure B9
pre_control %>% 
  filter(study=="Study 2" | (study=="Study 3" & treatment==0)) %>% 
  filter(!is.na(multi_item)) %>% 
  group_by(study, n, multi_item) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 3", 961, n), # adjust n in S3 to control group
         prop = freq/n) %>%
  ggplot(aes(x = multi_item, y = prop)) + 
  geom_bar(stat = "identity", width = 0.12, alpha = 0.7) + 
  facet_grid(study ~ .) +
  theme_params + 
  labs(x = "Multi-Item Attitude Certainty (1-9)", y = "Proportion") +
  theme(strip.text = element_text(size = 7),
        axis.text.y = element_text(size = 7),
        axis.title.y = element_text(size = 7),
        axis.text.x = element_text(size = 7),
        axis.title.x = element_text(size = 7))

# Figure B10
pre_control %>% 
  filter((study=="Study 4" | study=="Study 5") & treatment == 0) %>% 
  filter(!is.na(pw2p_scale)) %>% 
  group_by(study, n, pw2p_scale) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 4", 910, 649),
         prop = freq/n) %>%
  ggplot(aes(x = pw2p_scale, y = prop)) + 
  geom_bar(stat = "identity", alpha = 0.7) + 
  facet_grid(study ~ .) +
  theme_params + 
  labs(x = "Multi-Item Attitude Extremity (1-7)", y = "Proportion") +
  theme(strip.text = element_text(size = 7),
        axis.text.y = element_text(size = 7),
        axis.title.y = element_text(size = 7),
        axis.text.x = element_text(size = 7),
        axis.title.x = element_text(size = 7))

pre_control %>% 
  filter((study=="Study 4" | study=="Study 5") & treatment == 0) %>% 
  filter(!is.na(pw2p_defend)) %>% 
  group_by(study, n, pw2p_defend) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 4", 906, 643),
         prop = freq/n) %>%
  ggplot(aes(x = pw2p_defend, y = prop)) + 
  geom_bar(stat = "identity", alpha = 0.7) + 
  facet_grid(study ~ .) +
  theme_params + 
  labs(x = "Multi-Item Attitude Defense (1-7)", y = "Proportion") +
  theme(strip.text = element_text(size = 7),
        axis.text.y = element_text(size = 7),
        axis.title.y = element_text(size = 7),
        axis.text.x = element_text(size = 7),
        axis.title.x = element_text(size = 7))

# Table B19

## Line 1: Attitude Strength
pre_control %>% 
  filter(study!="Study 3" | (study=="Study 3" & treatment==0)) %>% 
  filter(!is.na(strength)) %>% 
  group_by(study, n, strength) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 3", 959, n), # adjust n in S3 to control group
         prop = freq/n) %>% 
  filter(strength == 7)

## Line 2: Attitude Certainty
pre_control %>% 
  filter(study!="Study 3" | (study=="Study 3" & treatment==0)) %>% 
  filter(!is.na(certainty)) %>% 
  group_by(study, n, certainty) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 3", 916, n), # adjust n in S3 to control group
         prop = freq/n) %>% 
  filter(certainty == 100)

## Line 3: Multi-Item Certainty
pre_control %>% 
  filter(study=="Study 2" | (study=="Study 3" & treatment==0)) %>% 
  filter(!is.na(multi_item)) %>% 
  group_by(study, n, multi_item) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 3", 961, n), # adjust n in S3 to control group
         prop = freq/n) %>% 
  filter(multi_item == 9)

## Line 4: Attitude Extremity
pre_control %>% 
  filter((study=="Study 4" | study=="Study 5") & treatment == 0) %>% 
  filter(!is.na(pw2p_scale)) %>% 
  group_by(study, n, pw2p_scale) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 4", 910, 649),
         prop = freq/n) %>% 
  filter(pw2p_scale == 7)

## Line 5: Attitude Defense
pre_control %>% 
  filter((study=="Study 4" | study=="Study 5") & treatment == 0) %>% 
  filter(!is.na(pw2p_defend)) %>% 
  group_by(study, n, pw2p_defend) %>%
  summarize(freq = n()) %>% 
  mutate(n = if_else(study == "Study 4", 906, 643),
         prop = freq/n) %>% 
  filter(pw2p_defend == 7)
```

## Appendix B.4 (Measurement Study)

```{r}
library(psych)
library(mirt)
library(qualtRics)
library(hrbrthemes)

df_fixed <- rio::import("data/measurement_study.csv")

# Reverse-coding specific items
reverse_coded_items <- c('pw2p_scale_8', 'pw2p_scale_9', 'pw2p_scale_10', 'pw2p_scale_11', 'pw2p_scale_12', 'pw2p_scale_13', 'pw2p_scale_14', 'pw2p_scale_15')
df_fixed[, reverse_coded_items] <- 8 - df_fixed[, reverse_coded_items]

# Calculate pw2p_scale mean
cols_to_average <- c('pw2p_scale_1', 'pw2p_scale_2', 'pw2p_scale_3', 'pw2p_scale_4', 'pw2p_scale_5', 'pw2p_scale_6', 'pw2p_scale_7', 'pw2p_scale_16', reverse_coded_items)

selected_items <- c("pw2p_scale_7", "pw2p_scale_14", "pw2p_scale_8", "pw2p_scale_1", 
                    "pw2p_scale_2", "pw2p_scale_9", "pw2p_scale_3", "pw2p_scale_6")

df_fixed$pw2p_scale_k <- rowMeans(df_fixed[, cols_to_average], na.rm = TRUE)
df_fixed$pw2p_score <- as.vector(fscores(mirt(df_fixed[, cols_to_average], 1, verbose = FALSE, itemtype = "gpcm")))
df_fixed$pw2p_scale_abb <- rowMeans(df_fixed[,selected_items], na.rm = T)

# Compute alpha for defense scale
df_fixed %>% 
 dplyr::select(contains("defend")) %>%
 psych::alpha(check.keys = TRUE) -> defense_scale

df_fixed$defense_score <- 2.60-as.vector(fscores(mirt(df_fixed %>% 
                                                       dplyr::select(contains("defend")), 1, verbose = FALSE,
                                                      itemtype = "gpcm")))

# Full scale
df_fixed %>%
 dplyr::select(contains("pw2p")) %>%
 select(-pw2p_scale_k, -pw2p_score, -pw2p_scale_abb) %>%
 psych::alpha(check.keys = T) -> full_scale

# Compute allocation differential
df_fixed <- df_fixed %>% 
 mutate(allocation_differential = post_allocate_1 - (post_allocate_2+post_allocate_3+post_allocate_4)/3)

# Update defense_scale scores
df_fixed$defense_scale <- 7 - defense_scale$scores

df_fixed$full_scale <- full_scale$scores

# Correlation analysis
cor_matrix <- cor(df_fixed[, c('defense_scale', 
                               'allocation_differential', 
                               'post_outcome_1', 'post_certainty_19',
                               'pw2p_scale_k', 'm_select')], use = "pairwise.complete.obs", 
                  method = "spearman")

print(cor_matrix)

# Label mapping
labels <- c('defense_scale' = 'Attitude\nDefense', 
            'allocation_differential' = 'Allocation\nDifferential', 
            'pw2p_scale_k' = 'Attitude\nExtremity', 
            'post_outcome_1' = 'Likert',
            'post_certainty_19' = 'Certainty',
            'm_select' = 'Conjoint\nSelection')

# Create the plot
library(ggcorrplot)

ggcorrplot(cor_matrix, hc.order = TRUE, type = "lower",
           lab = TRUE, lab_size = 3, method = "square", 
           colors = c("#6D9EC1", "white", "#E46726"), 
           outline.col = "grey") +
 ggtitle("Correlation Matrix") +
 scale_x_discrete(labels = labels) +
 scale_y_discrete(labels = labels) +
 theme_ipsum_rc(grid="XY") +
 theme(plot.title = element_text(hjust = 0.5, size=14, face="bold"),
       panel.grid.major = element_blank()) +
 labs(x = "", y = "")

# Min-Max scaling
df_fixed_minmax <- df_fixed %>% 
 mutate(across(c(defense_scale, allocation_differential, pw2p_scale_k, 
                 post_outcome_1, post_certainty_19,
                 full_scale), ~ ( . - min(.)) / (max(.) - min(.))))

# Calculate coefficient of variation for each variable
vars_to_check <- c("defense_scale", "allocation_differential", "pw2p_scale_k", "post_outcome_1", "post_certainty_19")

ceiling_effects_metrics <- function(dataframe) {
 library(dplyr)
 library(moments)
 
 skew_values <- apply(dataframe, 2, skewness, na.rm = TRUE)
 mean_values <- apply(dataframe, 2, mean, na.rm = TRUE)
 kurt_values <- apply(dataframe, 2, kurtosis, na.rm = TRUE)
 pct_max <- apply(dataframe, 2, function(x) sum(x == 1, na.rm = TRUE) / length(na.omit(x)) * 100)
 
 metrics_table <- tibble(
  Variable = names(dataframe),
  Skewness = round(skew_values,2),
  Kurtosis = round(kurt_values,2),
  Mean = round(mean_values,2),
  PctMax = round(pct_max,2)
 )
 
 return(metrics_table)
}

ceiling_metrics <- ceiling_effects_metrics(df_fixed_minmax[,vars_to_check])
ceiling_metrics %>% kableExtra::kable(format = "latex")

# Data in long format for plotting
df_long_minmax <- df_fixed_minmax %>% 
 select(defense_scale, allocation_differential, pw2p_scale_k, 
         m_select) %>% 
 gather(key = "variable", value = "value", -m_select)

# Create plot
ggplot(df_long_minmax %>%
       mutate(Variable = plyr::mapvalues(variable,
                                         c("allocation_differential",
                                           "defense_scale",
                                           "pw2p_scale_k"),
                                         c("Allocation Differential",
                                           "Attitude Defense",
                                           "Attitude Extremity"))), aes(x = value, y = m_select, color = Variable, group = variable)) +
 stat_smooth(method = "glm", 
             method.args = list(family = binomial(link = "logit")), se = FALSE) +
 labs(title = "Predicting Conjoint Selection", y = "Selection Rate") +
 theme_ipsum_rc() +
 theme(panel.grid.major = element_blank())

models <- list(
 lm(m_select ~ pw2p_scale_k, data = df_fixed_minmax),
 lm(m_select ~ defense_scale, data = df_fixed_minmax),
 lm(m_select ~ allocation_differential, data = df_fixed_minmax),
 lm(m_select ~ pw2p_scale_k + defense_scale + allocation_differential, data = df_fixed_minmax)
)

library(MASS)

modelsummary(models)

# Character count analysis
df_fixed_minmax$count <- stringr::str_count(df_fixed_minmax$introduction)

models <- list(glm.nb(count ~ defense_scale, data = df_fixed_minmax),
               glm.nb(count ~ allocation_differential, data = df_fixed_minmax), 
               glm.nb(count ~ pw2p_scale_k, data = df_fixed_minmax),
               glm.nb(count ~ pw2p_scale_k + defense_scale + allocation_differential, data = df_fixed_minmax))

modelsummary(models)

```

## Appendix B.5 (Lower Versus Higher Valence Counterarguments)

```{r}
source("scripts/clean_study5.R")

library(magrittr)
library(tidyverse)
library(stm)
library(hrbrthemes)
library(tidytext)

tf_idf_distinctive <- study5 %>%
 slice(3:1984) %>%
 transmute(`Lower Valence\nCounterargument` = softer_counterargument, 
           `Higher Valence\nCounterargument` = counterargument) %>%
 pivot_longer(cols = everything()) %>%
 unnest_tokens(word, value) %>%
 anti_join(stop_words) %>%
 count(name, word) %>%
 bind_tf_idf(word, name, n) %>%
 arrange(desc(tf_idf))

# Select top distinctive words for each category
top_words <- tf_idf_distinctive %>%
 group_by(name) %>%
 top_n(30, tf_idf) %>%
 ungroup() %>%
 mutate(word = reorder_within(word, tf_idf, name))

# Lollipop plot
ggplot(top_words, aes(x = tf_idf, y = word, color = name)) +
 geom_segment(aes(x = 0, xend = tf_idf, y = word, yend = word), size = 1) +
 geom_point(size = 4) +
 scale_y_reordered() +
 labs(title = "Top Distinctive Words by TF-IDF",
      x = "TF-IDF Score",
      y = "Word") +
 facet_wrap(. ~ name, scales = "free") +
 theme_ipsum_rc() +
 theme(axis.text.x = element_blank(),
       legend.position = "none") 

```


